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Abstract 


c 


The  effect  of  a  surface  nuclear  blast  on  a  deep  und'  r- 
ground  tunnel  was  studied  using  elastic-viscoplastic  theor  . 
Quasi-static  and  quasi-transient  loadings  were  analyzed  us  ng 
a  finite  element  model  developed  herein. 

The  computer  program  used  is  a  two-dimensional  plan' 
stress/plane  strain  code  using  constar. l  strain  triangle 
elements.  Constitutive  equations  proposed  by  Drucker  and 
Prager  were  incorporated  to  account  for  the  compressibi ii c 
inherent  to  rock  material.  Plasticity  was  generated  using 
the  Drucker-Prager  yield  criteria  incorporated  in  the  Maiv  rn 
flow  law.  Both  viscoplastic  and  elasto-plastic  material 
behavior  were  modeled  using  the  viscoplastic  theory. 

The  quasi-static  analysis  compared  excellently  to 
material  behavior  displayed  in  field  studies  and  experimun .el 
work.  Results  obtained  through  the  quasi-transient  analys 
give  an  appreciation  of  possible  plastic  zone  growth. 
Additionally,  a  modification  to  the  Drucker-Prager  equatio: 
was  proposed  and  found  to  show  less  plastic  strain  than  u: 
full  equations.  This  is  promising  when  considered  in  1  ig:.  o 
the  fact  that  the  Drucker-Prager  equations  have  been  known  • o 
allow  too  much  plasticity  at  high  pressure  levels. 
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VISCOPLASTIC  ANALYSIS  OF  A  CONTINUOUS 
CYLINDRICAL  OPENING  SURROUNDED 
BY  VOLCANIC  TUFF 


I.  Introduction 


Background 

Until  the  mid  1970s,  the  United  States  based  its  in  cr 
continental  ballistic  missiles  (ICBMs)  in  shallow  underyro  nd 
concrete  silos.  This  system  was  considered  sufficiently 
protective  to  allow  the  United  States  to  survive  a  "first, 
strike,"  thus,  giving  a  "second  strike"  capability. 

With  the  advent  of  higher  yield  warheads  couplea  vs  h 
increasingly  accurate  delivery  systems,  the  protection  pic 
vided  by  the  silo  system  became  suspect,  and  the  searci.  .o  a 
new  ICBM  basing  system  began.  One  of  the  ideas  which  cnnu 
out  of  the  search  is  the  Deep  Missile  Basing  System,  which 
envisions  a  network  of  tunnels  dispersed  throughout  a  J an. 
mass  called  a  mesa.  These  tunnels  are  buried  within  the 
mesa,  and  are  subjected  to  the  natural  in  situ  stresses, 
well  as  the  dynamic  stresses  resulting  from  the  spher:cai . 
propagating  stress  wave  generated  by  a  contact  or  near-.^.  i 
blast  (see  Fig  1)  (Ref  1). 

The  survivability  of  deep  undez'ground  tunnels  i. 


function  of  both  the  strength  of  the  employed  tunnel 
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Cutaway  of  Deep  Underground  Missile  Basing  Sysu 


hardening  measures  (concrete  liters,  steel  liners,  rock 
bolts,  etc.),  and  the  strength  of  the  surrounding  rock. 
Economic  tunnel  design  dictates  minimum  use  of  hardening 
techniques,  taking  advantage  instead  of  the  protective  roc 
surrounding  a  tunnel  (Ref  1).  Therefore,  the  state  of  str  = 
in  the  tunnel/rock  system  must  be  estimated. 

Approach 

Experimental  results  reveal  that  loading  near  tunne 
levels  from  a  surface  nuclear  blast  causes  failure  in  uni i  .c 
tunnels,  characterized  by  a  distressed  or  loosened  zone  of 
rock  at  the  springlines  (the  springline  is  perpendicular  t 
the  maximum  load  direction),  and  tension  cracks  at  the  croon 
Since  failure  normally  indicates  plastic  action,  an  elasti  - 
plastic  analysis  is  in  order. 

The  Drucker-Prager  plastic  mate. a  < 1  model  represent 
one  of  the  idealizations  of  real  rock  which  entails  apprcs 
elastic  constants,  a  yield  function  and  a  flow  rule.  'ihi\- 
important  characteristics  of  the  real  stress-strain  curve 
are  reflected  in  this  model.  First,  elastic  response  is 
pronounced  at  lower  loads.  Second,  as  load  rs  increaseu 
toward  ultimate,  the  actual  stress-strain  curve  has  a  Li 
bent  over  considerably  so  that  the  tangent  modulus  at  in.?, 
stage  is  a  small  fraction  of  the  initial  elastic  modules. 

The  perfectly  plastic  model  represent:  u  zero  tangent  n 


Finally,  the  plastic  behavior  of  rock  is  observed  oy  n^vi:. 
residual  strain  when  a  complete  unroadmg  takes  place  Ley 


the  elastic  range.  This  is  in  contrast  to  the  nonlinear 
elastic  idealization  where  unloading  follows  the  initial 
path  and  the  strain  is  fully  recoverable.  This  last  charac¬ 
teristic  gives  a  distinction  between  a  plastic  and  an  elas* ic 
rock  and,  in  the  most  fundamental  sense,  rock  has  been  knc.  n 
to  be  a  plastic  rather  than  an  elastic  material  (Ref  2) . 

The  rock/tunnel  system  will  be  modeled  using  the  fii.ite 
element  method,  with  an  in-house  code  named  VISCO.  The  ro>  k 
will  be  treated  for  simplicity  as  a  continuous  and  isotrop  c 
material.  For  the  depths  and  yields  of  interest  to  deep 
missile  basing  systems,  the  pressure  wave  from  a  nuclear  b  ast 
will  probably  have  a  significant  rise  time,  on  the  order  o 
1  to  10  transit  times,  and  be  of  fairly  long  duration,  100  ro 
200  transit  times;  in  effect,  making  a  quasi-static  (stead 
state)  approach  appropriate.  A  quasi-transient  analysis  v,  1 
also  be  done  using  appropriate  rise  times  for  each  load  lecl. 


Objective 

The  objective  of  this  thesis  is  to  analyze  the  strc 
around  the  tunnel  under  the  conditions  of  quasi-static  and 
quasi-transient  loading.  This  will  be  accomplished  using 
viscoplastic  theory  incorporating  the  Drucker-Prager  niatci  a 1 
model  developed  in  Chapter  II. 


II.  Theory 


Plasticity 

The  behavior  of  materials  beyond  the  elastic  range  .j 
termed  plastic,  and  can  be  characterized  by  an  incomplete 
instantaneous  recovery  of  strain  upon  removal  of  stress,  a 
change  in  strain  with  time  under  constant  load,  a  change  i  . 
stress  with  time  under  constant  deformation,  a  variation  o 
mechanical  properties  with  applied  stress,  or  a  variation  f 
mechanical  properties  with  temperature,  direction,  and  pcs  tion 
in  the  body.  Thus,  the  implication  is  that  the  theory  of 
elasticity  can  be  used  to  make  the  first  approximation  of 
stresses,  strains,  and  deformations  in  a  structure  under 
specified  loading,  and  that  various  inelastic  (plastic) 
theories  can  be  used  to  estimate  the  deviation  from  elast.i 
behavior  that  may  occur  and  to  estimate  the  total  loads  tr.  f  a 
given  structure  can  support  before  excessive  deformation, 
fracture,  or  disintegration  occurs  (Ref  4). 

Plasticity  in  Metals 

Based  on  the  above  paragraph,  total  strain  in  a  mat  : ial 
can  be  written  as 


total 


elastic 


plastic 


Lij 


+  c  .  . 
ID 


V.  ) 


The  onset  of  plastic  deformation  (or  yielding)  is  govua> 
a  yield  criterion,  and  post  yield  deformation  occurs  ... 


greatly  reduced  material  stiffness.  A  yield  criterion  o,. 
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where  f  is  some  rune _J on 


written  in  the  form  f(a^j)  =  ^  • 

and  k  is  a  material  parameter  to  be  determined  experimen 
Experimental  observations  by  Bridgman  (F<ef  3)  indicate  tlia 
plastic  deformation  of  metals  is  essentially  independent  o 
hydrostatic  pressure.  Therefore,  yield  criteria  can  now  b 
written  as  3  k  ,  where  and  are  the  se 

and  third  invariants  of  the  deviatoric  stresses 
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ij 


1 
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One  of  the  simplest  yield  criterion  is  the  Von  Misc 
T 

criterion  ( J ')  =  k  shown  in  Fig  2,  which  depicts  the  vi 

surface  and  the  locus  of  points  upon  which  the  principal  n 
stresses  sum  to  zero  (n  plane) .  This  may  also  be  written 

i_ 

a  =  JTtJ')  2  =  J7k  U 

where 


and  o  is  termed  the  effective  stress,  a  stress  invariant 

measures  the  incremental  change  of  the  yield  surface  to 

overstress.  For  the  case  of  uniaxial  tension  |on  - 

]3k  is  the  uniaxial  yield  stress.  For  yielding  to  occur, 
1  _L 

( J 2 )  must  reach  the  value  of  (Jp  ir;  a  uniaxial  \ 

test.  'J'he  Von  Mises  yield  criterion  implies  that,  il 
of  stress  lies  on  the  yield  surface,  plastic  riow  occur.,, 
the  state  of  stress  is  within  the  yieiu  surface,  elastic  . 


ally. 

end 

■-) 

Ld 
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•  i 
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>) 

which 


a  ins 


occur.  During  plastic  straining,  this  yield  surface  expands 
normal  to  itself  (normality  condition) ,  thus,  enforcing 
isotropic  hardening  within  the  plastic  zone. 


After  the  onset  of  plasticity,  a  so  called  equivalent 
yield  stress  function  can  be  evaluated  such  that  a  gauge  on 
plastic  flow  can  be  observed,  and  which  increases  as  a  function 
of  the  plastic- strain.  The  equivalent  yield  stress  can  be 
written  as 

Y  =  oy  +  H'cp  (2.5) 

where  is  the  uniaxial  yield  stress,  H"  is  the  slope  of 

the  strain-hardening  portion  of  stress-strain  curve,  and  ep 
is  the  plastic  strain  equivalent  (Ref  o)  .  Y  is  used  as  a 
universal  stress-strain  curve  that  governs  the  material's 
uniaxial  inelastic  behavior  during  loading,  and  is  valid  only 
when  used  with  monotonically  increasing  loads.  Once  the  onset 
of  plasticity  has  occurred,  increased  plasticity  will  result 
only  if  o  exceeds  Y 

Plasticity  in  Rock 

Rock  differs  from  metal  in  terms  of  plasticity  because 
rock  is  a  compressible  material.  While  metals  have  one  stress- 
strain  curve,  rocks  have  several  which  are  a  function  of  the 
confining  stress  ~  °3  •  An  increase  in  confining  pres¬ 

sure  significantly  increases  the  compressive  strength  of  rock. 
Most  ductile  materials  have  a  fairly  long  plastic  range  before 

fracture  occurs,  where  as  brittle  materials  such  as  most 
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rocks  have  a  very  limited  plastic  range,  if  any,  for  uniaxial 
tension  or  compression.  However,  for  triaxial  loading,  most 
brittle  materials  tend  to  show  a  plastic  region  where  permanent 
deformation  occurs.  A  yield  criterion  which  accounts  for 
hydrostatic  stresses  is  the  Drucker-Prager  yield  criterion 
(Ref  2) ,  which  modifies  Von  Mises  by  including  ,  the  first 

stress  invariant  (J^  =  +  O22  +  O33)  •  Following  the  same 

analogy  as  with  the  Von  Mises  (including  the  normality  condi¬ 
tion)  ,  the  Drucker-Prager  criterion  car  be  written  as 

1 

aJj  +  (J2'  )  =  k"  (2.6) 

Therefore,  plastic  flow  begins  when  the  effective  stress 
o  =  aJ^  +  (j')  2  exceeds  k"  ;  k"  can  be  thought  of  as  the 
yield  stress,  which  in  this  case  is  a  constant.  The  constant 
a  describes  the  dilatation  rate  of  soil,  dilatant  behavior 
being  a  volume  increase  under  compression  (Refs  6,  7,  17). 
Digressing  somewhat,  the  Drucker-Prager  yield  criterion  is  an 
approximation  to  another  criterion,  the  Mohr-Coulomb  yield 
criterion.  Note  in  Fig  3  that  two  Drucker-Prager  yield 
surfaces  are  required  to  satisfy  Mohr-Coulomb  at  the  apices. 
Varying  the  radii  is  accomplished  by  using  two  equations  for 
a  and  k'  .  For  the  outer  apices  (BDF) 


2  sinQ _ 

J  3  (3-sin<|>) 


and  k' 


6c  (cosl) 


J3  (3-sinl) 


(2.7) 


and  for  the  inner  apices  (ACE) 


ot 


and  k 


6c  (cos<{>) 

(3  +  sin4>) 


(2.8) 


I 52 


* 

k  •**% 

B  (•' 


2  sinj> _ 

JT  (3  +  sin<j>) 


C  is  the  cohesion  coefficient,  a  measure  of  the  shear  resistance 
between  any  two  particles  in  a  rock  mass;  shear  resistance 
is  the  combined  effect  of  cohesive  resistance  and  frictional 
resistance.  <p  is  the  angle  of  internal  friction,  and  is 
related  to  the ‘ fracture  plane  of  a  specimen  occurring  under 
a  triaxial  state  of  stress.  Studies  dealing  with  softer 
rocks  and  soils  (Refs  10,  11)  use  Eq  (2.8),  while  Eq  (2.7) 
apply  to  stiffer  materials  such  as  harder  rocks  or  concrete. 


Viscoplasticity 

A  time-oriented  viscoplastic  finite  element  code  called 
VISCO  is  used  in  this  study.  Time  rate  effects  are  always 
present  to  some  degree  in  all  inelastic  deformations.  Whether 
or  not  their  inclusion  has  a  significant  influence  on  the 
prediction  of  the  material  behavior  depends  upon  several 
factors.  In  the  study  of  structural  components  under  static 
loading  at  normal  temperatures,  it  is  accepted  that  time  rate 
effects  are  generally  not  important  and  the  conventional 
theory  of  plasticity,  as  described  above,  models  the  behavior 
adequately.  However,  some  materials,  especially  at  high 
temperatures,  exhibit  simultaneously  the  phenomena  of  creep 
and  viscoplasticity.  Creep  can  be  defined  as  a  redistribution 
uf  stress  and/or  strain  with  time  under  elastic  material 
response  while  viscoplasticity  is  a  time  dependent  plastic 
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deformation  (Ref  5) .  A  viscoplastic  analysis  increments  time 
under  full  loading,  while  an  elastic  plastic  analysis  incre¬ 
ments  only  the  load. 

The  concept  of  viscoplastic  material  behavior  is  detailed 
through  the  use  of  the  one-dimensional  rheological  model  in 
Fig  4.  The  friction  slider  develops  a  stress  Op  , 
becoming  active  only  if  a  >  Y  ,  where  Y  is  as  above.  The 
excess  stress  =  a  -  a^  is  carried  by  the  viscous  dashpot . 
Instantaneous  elastic  response  is  carried  via  the  spring. 

This  model,  specifically  the  dashpot,  allows  the  stress  level 
to  instantaneously  exceed  the  value  prescribed  by  plasticity 
theory,  the  solution  tending  to  a  steady  state  condition  with 
time  (Ref  5) . 

Viscoplastic  strain  in  rate  form  can  be  written  as 


VP 

ij =  y 

<  <f(F)  > 

3  F 

3a 

(2.9) 

<t (F) > 

=  <MF) 

for 

F  >  0 

(2.10) 

<$ (F) > 

=  0 

for 

F  <  0 

(2.11) 

>  is  a  fluidity  parameter  controlling  ^Listic  flow  rate,  and 

0  F 

is  normally  derived  from  uniaxial  yield  tests,  is  the  flow 

vector  a  for  the  associated  plasticity  which  governs  strain 
in  the  plastic  region.  The  form  of  : (F)  used  in  VISCO  is 


where  F  corresponds  to  o  (effective  stress)  and  Fq 
corresponds  to  Y  (equivalent  yield  stress)  as  described 
previously.  Using  the  Von  Mises  flow  vector  (Ref  5)  results 
in  the  VISCO  definition  of  the  viscoplastic  strain  rate 


;2-rr2-il2l£ 

3  L*  J  2  pr 


(2.13) 


where  are  the  deviatoric  stresses. 

Using  the  Drucker-Prager  flow  vector  and  recalling  that 

i 

Y  =  k^  and  cr  =  aJ^  +  ( J ')  2  the  viscoplastic  strain  rate 
can  be  written  as  (Ref  6) 


5  ■  <  E-  ■]  [•„] 


(2.14) 


where 


w  .  . 

.  .  =  a  +  — ^ — -_i_  when  i  =  j 


2(J2'H 


(2.15) 


a  _  - j 

1]  (J2)? 


when  i  /  j 


(2.16) 


The  generated  plastic  strain  rates  are  multiplied  by  a  time 
increment  to  yield  plastic  strains.  As  noted  earlier,  a 
steady  state  analysis  will  also  be  done.  Under  this  condition 
the  key  to  using  VISCO  is  the  fluidity  parameter  y 
Perzyna  showed  that  the  constitutive  equations  of  plastic 
flow  theory  result  as  a  limiting  case  from  the  more  general 


constitutive  equations  of  the  viscoplastic  theory  as '  y 
approaches  infinity  (Ref  6) .  The  definition  of  steady  state 
considering  a  finite  y  in  conjunction  with  Eq  (2.14)  will 
be  discussed  later. 


III.  VISCO  Computer  Program 


VISCO  is  a  plane  strain/plane  stress  finite  element 
code  utilizing  constant  strain  triangle  elements.  Because  the 
program  was  initially  designed  for  crack  growth  analysis,  it 
uses  a  Gauss-Seidel  iterative  solution  technique  employing 
the  initial  stress  method  (Ref  10)  to  solve  the  equilibrium 
equation 


[k] { d )  =  {R}  (3.1) 

where  [k]  is  normally  highly  banded  and  stored  in  compacted 
form  (retaining  only  non-zero  terms) ,  and  {d}  is  the  nodal 
displacement  vector. 

As  indicated  above,  once  an  element  yields,  the  equiva¬ 
lent  yield  stress  increases  as  a  function  of  the  plastic 
strain.  This  can  be  thought  of  as  nonlinearity  induced  through 
material  property  changes.  Two  procedures  exist  to  account  for 
the  effects  of  plasticity:  one  procedure  reformulates  the 
matrix  each  time  an  element  experiences  plastic  deformation 
and  is,  therefore,  both  time  consuming  and  expensive.  The 
second  procedure  is  called  the  residual  force  method  in  which 
the  elastic  stiffness  matrix  is  used  throughout  the  analysis 
and  the  effects  of  plasticity  are  treated  as  an  applied  load 


to  be  added  to  {R; 

If  the  elasto-plastic  stiffness  matrix  is  written  as 


► 


r 


H  <• 

■ 


♦ 


where  [kJC  is  the  correction  to  the  elastic  stiffness  matrix 
at  a  given  stage  in  the  analysis.  Equation  (3.1)  for  a 
particular  load  increment  i  becomes 

<[k]E  -  [ki]C){Ad.}  =  {AR.  }  (3.3) 

and 

[k]E  {Adi>  =  {ARi)  +  [ki]C  {Ad..}  (3.4) 

Defining 

{Ad;.}  =  {Adi)E  +  ^di>P  (3.5) 

where  {Ad^}  is  the  incremental  displacement 
then 

[k]E  {{Ad  }E  +  {Ad.  f)  »  {AR. }  +  [k,]C  ({Ad.  }E 
J  1  i  l  1  x 

+  {Adi)P)  (3.6) 

from  the  elastic  analysis 

[k]E  { Ad± }E  =  (ARi)  (3.7) 

so  that 

[k]E  { Ad^  } P  =  [ k i]C  (Ud.}E  +  (AdilP)-AP  (3.8) 

In  short,  the  additional  displacements  due  to  plasticity  can 
be  solved  by  applying  a  corrective  fo*  ;e  AP  to  the  elastic 
system  (Ref  11). 
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The  residual  force  method  can  be  utilized  even  for  a 
viscoplastic  solution  by  incorporating  a  time  step  algorithm 
as  follows  (Ref  9) : 

1.  Increment  time  so  that  t  =  t^-^  +  dt1 

2.  Compute  the  plastic  strain  increment  as  described  in 

above  Eqs  (2.14),  (2.15),  =  { e  ^  } 1  dt^  and  add  to  the 

preceding  plastic  strain. 

3.  Compute  the  plastic  load  vector 

{p}1'1  =  /  [b]  [k]  (h.p.  )1  dvol 

vol  3 

where  [b]  is  the  constant  strain  triangle  strain-displace¬ 
ment  matrix. 

4.  Compute  the  external  load  vector 

0 

{R}1  =  {R}1"1  +  (R}dt£ 

5.  Compute  nodal  displacements 

(d}1  =  [k]"1  ({R}1  +  {P  }i_1) 

6.  Compute  current  total  strain 

{eij}i  =  [B]  {d}1 

7.  Compute  current  stress 

tcij,i  “  W  (u^)1  - 
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IV.  Rock  Mechanics 


n 


I 


I 


9  o 


fc 


The  state  of  stress  at  a  point  strongly  influences 
the  strength,  stiffness,  ductility  and  creep  properties  of 
rock.  Figure  5  shows  variation  in  confining  pressure, 

versus  axial  strain  for  Sandston  (Ref  13)  .  Here,  as 
in  the  remainder  of  this  report  (unless  otherwise  noted) ,  all 
displayed  stresses  are  compressive.  The  increase  in  confining 
pressure  results  in  a  significant  increase  in  ultimate  strength 
as  well  as  a  larger  strain  at  failure.  Also,  plastic  strains 
before  ultimate  failure  (fracture)  increase  with  confining 
pressure.  The  Drucker-Prager  yield  criterion  accounts  for 
confining  pressure  by  including  the  first  stress  invariant 
(J^  in  Eq  (2.6)),  and  also  through  the  two  material  parameters 
C  and  4)  defined  above.  C  and  are  normally  derived 

from  triaxial  test  data  using  Mohr's  circles  as  shown  in 
Fig  6.  Note  that  ck  values  are  stresses  at  failure;  in 
essence,  the  circles  represent  the  material  failure  surface  as 
a  function  of  confining  pressure.  An  alternative  method  was 
used  here,  in  which  C  was  derived  from  a  known  $  value 
and  uniaxial  test  data.  Using  the  spherical  stress  component 


x  - 


then 


I  =  +  JY  +-.°_  and  recalling  that  a 

3  XY  * 

Ov  Ox 

for  a  uniaxial  test  I  =  _  and  ci  ,  =  which  yields  a 

3  xy  ^ 

point  on  a  plot  of  o  versus  I  (not  shown) .  Using  a 

xy 

sponsor  supplied  (from  an  unpublished  report)  slope  of  ;  -  18 

a  o  intercept  is  determined  and  is  used  as  the  value  lor 
xy 


o 
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V.  Modelirv 


Mesh  Size 

Since  the  size  of  the  modeled  section  for  this  study 
had  the  potential  to  become  very  large,  it  was  necessary  to 
discover  early  a  satisfactory  mesh  size  limit.  Mesh  grading 
was  also  considered  in  conjunction  with  appropriate  aspect 
ratios,  where  aspect  ratio  (AR)  is  the  ratio  of  triangle  height 
to  width  (not  greater  than  AR=3)  (Ref  14) . 

Timoshenko  (Ref  15)  gives  a  closed  form  elastic  solu¬ 
tion  to  the  problem  of  a  hole  in  a  plate  under  uniaxial  tension 
or  compression,  based  on  a  ratio  of  d  to  a  ,  where  d  is 
the  plate  half-width  and  a  is  the  radius  of  the  hole  (Pig  7)  . 
Values  of  the  stress  concentration  factor  kc  at  the  crown 
and  springline  were  used  as  the  basis  for  determining  the 
fineness  of  the  mesh  around  the  cutout.  The  radius  of  the 
hole  is  problem  dependent;  thus,  determining  d  for  any  given 
d/a  ratio.  Arbitrarily  selecting  a  d/a  ratio  of  5,  with 
a  given  as  60  inches,  yields  a  d  value  of  300  inches  for 
the  half-plate  width.  Using  an  element  size  of  6  inches  by 
12  inches  (AR=2)  near  the  hole  yielded  values  of  2.97  for 
at  the  springline,  and  -1.01  for  at  the  crown.  Comparing 

to  Timoshenko's  values  of  3.08  for  and  -1.0  for  ox  , 

gave  acceptable  errors  of  3.56%  and  1.0%  respectively. 


Model  Size 

The  in  situ  state  of  stress  in  a  geological  mass  is  the 
natural  state  of  stress  at  any  given  point.  Baker,  et  al 
(Ref  11)  used  a  finite  element  model  employing  both  elastic 
and  elasto-plastic  theory  to  analyze  a  deep  tunnel  in  a  rock 
mass.  Setting  both  horizontal  and  vertical  boundaries  at 
10  radii  (600  inches)  gave  elastic  in  situ  stresses  in  agree¬ 
ment  with  Baker's (maximum  errors  of  4.0%  for  a  and  6%  for 
ox  on  a  horizontal  axis  through  the  tunnel,  see  Fig  8). 
Summarizing,  the  model  is  600  inches  square  with  323  nodes  and 
576  elements.  The  smallest  elements  arc  next  to  the  cutout  and 
are  6  inches  by  12  inches.  As  Fig  9  shows,  boundaries  not 
loaded  are  fixed  to  move  in  only  one  direction. 


26 


.00  40. CO  80.00  120.00 


o 

o 


o 


J  .  uu 

RflDI  I 


Fig  8.  ux  and  Oy  Comparison  Between  Published 

Results  (Ref  11)  and  576 
Element  Model 


27 


/  SJ9\ 

/  503\ 

V  “’\ 

/ 

411^ 

y 

395  ' 

\  /  359  \  / 

3?3  \ 

A^« 

\  974 

/\.  538 

/N.  SO? 

/'Nv  466  y 

\ 

430 

/\ 

394 

A\  is*  A\ 

i 

32’  / 

\ 

sn\ 

/  6”\ 

/  S0A 

/  46S\  / 

42  9N 

v 

TJ3N 

V 

/3Sl\  A?‘\ 

266 

29  1  A 

c  ^ - - "n 

\  51?  1 

5lV\  / 

\  536  1 
S3fi\  / 

r\  ^uo  / 
♦  9JV  / 

\  464  /\ 

483\  /  42V" 

428 

V 

A 

392 

V 

4X 

356  /  > \  370  / 

A'A/.( 

28S  y 

/ 

/  252 

970  / 

\  “A 

v  498  /\ 

462  /\  426 

7\ 

330 

/V 

3S4  / 

A  5,6 A 

A  j5° 

2Sl^"s 

y 

S6J\  /  SjK  /  48< 


S32  7\  <96  7\  <60  7\  424  7V  588  A  is,/\  316 


S6>\^ 

1 S3> 

y^y  4s\^y  4S?^ 

/ 42  \y  3fl  \/ 

36  )\  . 

N.  S66 

f\  S30  /\  4  04  /\  4S8 / 

/  5>V  /  4  a%  /  4S\  /4 

sit  r\  ,86/  \  3 

2l\  /3A  /  349\ 

50 y\ 
/113\ 

'\  S6  4 

56V  1 

r\6’e 

/XAAA7AA 

312/ 

A16 

i  y  ‘‘y, 

\W 

^7 

\  4J0/\4S4/\4l9/ 

'YYY 

\3<V\/\  34  6/1  31(3/ 
\/ J\/^U 

Tts/1 

/?4  0 

A 2,9/\ 218 

2  4  3  >A  2 1 «  JIS\ 
A  212  21  jV  y'  _ 


IT \  smV^vV  4iy  \  *  s7A7  ‘  Va3  *  A  3  4  tX 308/  7tT7 

/4^\^/4 s\^/4  ^  /&r  yfi 2  y'^ 

7t  7 ,  A  rj  \Av 6 1  A? 1  vAiAtAioV  7n?  A- 


A^A*08  ?09A  A 116 

/?06  ?#2^>Ai3  4 
2A>Al7?  \T)\ 


^  <A  Sk  •"''n  Nv\8  ^ 


j.tmx/v  V-V-AA  AiA'  >%>< 

jpm  fA  <a  AA3'1 
t  ^  <a  *  - 1 


K30  iK1 

•J3\ 

In  9s.  A* 

cTT 

cA  \6(i 
5K 

S7\ 

A* .  —  ~ 

'TTb  x2o 

\22 

1 7.  iK  | 

:’i\ 

^V!  .> 


Tv  Vs 
1 1  »\  !  A 

k - A - 


,'j\  I  71 


Fig  9.  Finite  Element  Model:  323  Node.;,  57b  Elements 
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One- Dimensional 


To  verify  the  changes  made  in  VISCO  to  employ  the  Drucker-Prager 
elastic-perfectly  plastic  material  model,  a  one-dimensional  viscoplastic 
analysis  was  done  using  a  plane  stress,  simple  tension,  two-element  model. 
Owen  (Ref  5)  notes  that  for  this  case,  a  stress  unbalance  will  always  exist 
and  steady  state  conditions  cannot  be  achieved  using  the  initial  strain 
approach  of  VISCO  as  would  be  expected.  In  equation  form 

e’vp  =  yO1  -  (oy  +  He^p)]  (6.1) 

but  since  H"  =  0 

4  ■  yCo1  -  °y]  (6.2) 

where  (a1  -  ay)  is  the  stress  unbalance.  Using  the  algorithm  described  in 
Chapter  III  where  stress  strain  relationships  are  in  terms  of  the  elastic 
stiffness,  A a  or  the  increment  of  stress  between  time  steps  can  be  written 
in  the  following  for'": 


Aa  =  E(Ae1'  -  Ae^) 


(6.3) 


The  two  element  case  was  run  over  a  period  of  ten  seconds  and  the  stress 
increment  approaches  zero  as  it  should  so  that  response  is  contributed  to 
the  time  interval  size.  For 


A a  =  E ( Ae '  -  Ae1 Vp)  =  0 


(6.4) 


Ae1'  -  Ae\ 


(6.5) 


Fig.  10  shows  a  delay  in  acquiring  the  zero  Aa  Ae  t  Aevp).  This  can  be 
attributed  to  the  size  of  the  time  interval  incorporated  into  the  viscoplastic 
equation.  The  smaller  the  time  step  the  closer  the  differential  strains. 

The  Visco  model  responded  as  indicated  by  Owen  in  that  there  is  always  a 
stress  unbalance,  flarcal  (Ref  8)  has  shown  that  a  time  independent  approach 
(non-visco)  would  produce  a  singularity  for  the  same  conditions. 


Two-Dimensional 


The  Baker,  et  al  (Ref  11)  elastic-plastic  tunnel  analysis  was  used 
as  another  verification  of  the  Drucker-Prager  model  used  viscoplastically. 
Previous  analysis  showed  that  a  stress  concentration  exists  around  the  cutout 
when  analyzed  elastically.  Note  in  Fig  11  the  comparison  between  the 
stresses  developed  in  the  elastic  solution  and  those  developed  in  the  elasto- 
plastic  solution  along  a  horizontal  section  through  the  tunnel.  Steady  state 
conditions  were  achieved  using  a  fluidity  parameter  y  of  1.0  in2/lbf-sec 
Steady  state  is  indicated  by  a  zero  strain  rate  (Eq  2.14)  in  all  elements. 

The  results  illustrate  the  fundamental  principle  that  plastic  flow  tends  to 
relieve  the  high  stress  concentrations  that  would  ordinarily  develop 
in  perfectly  elastic  materials.  Because  stresses  are  redistributed, 
note  that  the  viscoplastic  stresses  exceed  the  elastic  stresses 
as  distance  from  the  cutout  increases.  These 


3 


and 


results  compare  well  with  Baker's;  within  15%  for 
within  5%  for  ox 

Due  to  the  high  cost  and  lengthy  computer  time  required 
for  solution  of  the  equilibrium  equation  when  using  the  576 
element  model,  a  new  model  was  developed  (Fig  12)  containing 
only  238  elements.  The  576  element  model  (fine  mesh)  required 
computer  time  j(CPU  time)  of  700  seconds  at  a  cost  of  $52.34 
to  solve  Baker's  elastic-plastic  problem,  while  the  238  element 
model  (course  mesh)  required  285  seconds  at  a  cost  of  $20.38. 

A  comparison  between  the  two  meshes  in  terms  of  stresses 
and  a  along  three  representative  radials  is  detailed  in 
Figs  13-18;  note  that  the  two  meshes  do  not  vary  by  more  than 
10%.  The  accuracy  developed  by  the  coarser  mesh  versus  its 
lower  cost,  dictated  use  of  the  model  ior  the  remainder  of  the 
study. 

The  final  comparison  designed  to  validate  both  program 
and  model  contrasts  material  movement  generated  experimentally 
(Ref  1)  versus  material  movement  predicted  using  the  theory  and 
model  developed  thus  far.  Using  loads  that  approximate  loading 
at  tunnel  levels  from  a  nuclear  blast,  resulted  in  nodal 
displacements  presented  in  Fig  19.  Note  that  in  the  vertical 
direction,  most  material  moves  up  away  from  the  cutout  while  the 
node  at  the  springline  moves  in  the  opposite  direction  because 
the  cutout  provides  no  support  (all  cutout  nodes  collapse 
inward).  This  compares  excellently  with  experimental  results. 


COARSE 


Fine  mesh--576  element  model 
Coarse  mesh — 23d  element  model 

NOTd:  Stress  profile  alonj  a 

17°  ray  from  the  horizontal. 


Fine  mesh — 57 C  element  model 
Coarse  mesh  —  23d  elemeit  model 

NOTE:  Stress  profile  alonj 
4  5°  ray  from  the  horizontal. 
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VII.  Results  and  Discussion 


Quasi-Static  Analysis 

Two  different  loadings  are  considered  in  this  section 
and  in  the  following  section,  representing  either  varying 
nuclear  weapon  yields,  varying  tunnel  depths,  or  both.  In 
situ  loads  were  also  considered,  and  are  merely  added  to  the 
pressure  wave  loading.  For  convenience,  the  lower  magnitude 
loading  will  be  referred  to  as  1/4  kbar  (1  kbar  =  15,000  psi) , 
the  higher  loading  as  1/2  kbar  loading.  The  weapons  effect 
loading  shown  (Figs.  15,  16)  can  be  attributed  to  wave  propa¬ 
gation  studies  as  indicated  in  Reference  1. 

Figures  20  and  21  show  the  plastic  zone  developed  for 
1/4  and  1/2  kbar  loading  respectively,  using  a  y  value  of 
1.0.  The  reader  should  recognize  that  even  though  a  steady 
state  solution  is  being  carried  out,  a  relatively  small  value 
for  y  is  being  used.  This  value  simply  means  that  steady 
state  is  obtained  after  a  longer  period  of  time  compared  to 
larger  y's  .  Note  that  the  plastic  zone  envelopes  the  entire 
cutout  for  the  larger  load  and  that  the  vertical  applied 
stress  (3475  psi)  for  this  larger  load  has  increased  in  greater 
proportion  than  the  horizontal  stress  between  the  two  cases. 
Figures  22  and  23  display  tunnel  deformation  for  the  1/4  and 
1/2  kbar  cases  respectively.  It  can  be  seen  that  the  1/2  near 
case  experiences  much  larger  deformations  than  the  1/4  kbar 
case.  Also  note  that  the  major  deflection  occurs  at  the  crown 
for  the  1/4  kbar  loading,  and  at  the  springline  for  the  1/2  kbe 


3475  psi 


Fig  21.  Steady  State  Plastic  Zone--l/2  kbar 


case.  The  much  larger  deformations  for  the  1/2  kbar  loading 
are  due  to  both  the  increased  load  and  the  fact  that  the  entire 
tunnel  is  experiencing  plastic  deformation.  Figure  19  shows 
nodal  displacements  in  the  model  for  the  1/4  kbar  case;  dis¬ 
placements  are  similar  for  1/2  kbar  loading  (Fig  24).  Both 
show  springline  area  displacements  consistent  with  exper imenta- 
tion  (Ref  1) . 

Figure  25  compares  plastic  zone  stresses  along  a  45° 
radial  between  the  two  cases,  where  positive  values  denote 
ccmpression.  The  effective  stress  a  for  all  elements  in 
the  plastic  region  is  the  same  (976.4  psi) ,  despite  the  much 
larger  stresses  experienced  for  1/2  kbar  loading.  It  was 
noted  earlier  that  plastic  flow  redistributes  high  stress 
concentrations.  This  redistribution  stresses  is  clearly 
seen  here  because  the  much  larger  stresses  developed  in  the 
1/2  kbar  case  versus  the  1/4  kbar  case  still  produce  an 
effective  stress  equal  to  the  constant  yield  stress  k" 

Note  the  variance  in  nodal  displacements  displayed  in  Figures 
19  and  24.  Because  the  Drucker-Prager  yield  criterion  accounts 
tor  hydrostatic  stresses  (compressibility) ,  a  large  amount  of 
material  movement  through  compaction  and  crushing  occurs. 

Thus,  the  effective  stress  equals  the  yield  stress  for  two 
reasons:  relative  displacements  are  »  t  as  large  as  they 

initially  appear,  and  the  model  will  only  allow  the  affective 
stress  to  equal  the  yield  stress  in  the  plastic  zone  for  steady 
state  conditions  (contrasting  the  results  obtained  below) . 


However,  this  yield  stress  does  not  correspond  to  a  unique 
strain  level,  as  shown  in  Figures  26  and  27  which  display 
effective  plastic  strain  (Eq  3.2  below)  contours  in  the  plastic 
zone  for  the  1/4  kbar  case  and  the  1/2  kbar  case  respectively. 

The  larger  effective  strains  in  the  1/2  kbar  case  (at  least  6 
times  greater)  indicate  more  straining  occurring,  which  explains 
that  case's  larger  plastic  zone. 

Twenty-four  elements  in  the  1/2  kbar  loading  model 
experienced  tension  in  one  or  both  component  directions  (x  or 
y) .  These  elements  were  all  within  27  inches  of  the  cutout  and 
on  radials  from  11°  to  90°  from  the  horizontal  (note  that 
Fig  25  shows  no  tension  because  plotted  stresses  are  averaged 
between  two  adjacent  elements  and  compression  always  dominated 
tension).  However,  9  of  the  24  tension  elements  were  concentrated 
between  the  85°  and  90°  radials.  This  concentration  of  tensile 
stresses  is  appropriate  when  viewed  in  terms  of  material 
movement  at  the  springline  depicted  in  Fig  24.  The  1/4  kbar 
case  also  had  tension  elements,  although  fewer  and  confined  to 
radials  from  75°  to  85°.  Comparing  the  two  cases,  strain 

components,  as  well  as  effective  strains,  were  on  the  ordoi  of 
4  to  5  times  larger  for  tension  elements  in  the  1/2  kbar  cases 
versus  1/4  kbar.  Strains  this  large  could  be  indicative  ot 
tensile  material  failure. 

(juasi-Transient  Analysis 


In  order  to  mote  accurately  represent  the  physical  wave 
propagation  through  the  tunnel,  a  triangular  loading  function 


(as  indicated  in  Ref  1)  was  incorporated  into  the  analysis. 

Two  different  load  rates  (supplied  by  the  sponsor)  were  used, 
which  are  a  function  of  the  loading  levels. 

Figures  28-31  show  plastic  zone  development  for  quasi¬ 
transient  loading.  Note  the  differing  values  of  the  fluidity 
parameter  y  (0.1  and  1.0)  chosen  to  provide  a  lower  and  upper 
bound  for  possible  y  values  to  be  determined  experimentally. 
Recall  that  y  =  1.0  was  used  in  the  steady  state  analysis. 

The  same  y  value  is  now  being  used  to  study  a  transient 
state.  This  can  occur  due  to  the  unique  set  of  equations 
present  in  the  viscoplastic  solution.  As  employed  here,  > 
is  truly  acting  as  a  time  dependent  function,  rather  than  as 
a  so  called  dummy  variable  as  it  was  in  the  steady  state 
approach.  Each  figure  shows  the  load  history;  P  vertical 

(1000  psi  in  situ  load)  does  not  equal  P  .  horizontal 

min 

(467  psi  in  situ  load) .  Loads  shown  are  maximum  loads  (in 
situ  plus  weapons  effect) .  Comparison  with  the  static  results 
(Figs  20,  21)  reveal  significantly  smaller  plastic  zones  for 
both  load  levels,  because  the  short  times  modeled  do  not  allow 
as  much  plastic  flow  as  in  the  steady  state  analysis.  The 
size  of  the  1/4  kbar  transient  plastic  zone  was  27%  and  47  u 
of  the  static  plastic  zone  for  y  =  0.1  and  ,  =  1.0 
respectively.  For  the  1/2  kbar  case,  1  5'i  and  43%  of  the 
static  plastic  zone  for  y  =  0.1  and  ,  -  1.0  respectively. 
The  plastic  zones  for  =1.0  (Figs  29,  31)  are  larger  than 
the  zones  for  y  =  0.1  (Figs  28,  30)  because  )  =0.1  pro¬ 
vides  less  plastic  flow  as  dictated  by  Eq  (2.14). 
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Mote;  Shaded  area  is  plaR- 
tic  zone  at  t=.015  sec. 
Total  plastic  zone  enclosed 
by  dark  line. 


] 


Mote;  Shaded  area  is  plas¬ 
tic  zone  at  t=.01  sec. 

Total  plastic  zone  enclosed 
by  dark  line. 


Recall  that  viscoplastic  material  behavior  as  modeled 
by  che  one-dimensional  rheological  model  (Fig  4)  allows  the 
stress  level  to  instantaneously  exceed  che  yield  stress  due 
to  the  overstress  model  employed.  Effective  stresses  o' 
greater  than  the  yield  stress  are  shown  in  Figs  32  and  33, 
which  are  plots  of  o  versus  time  for  an  element  in  the  plastic 
zone.  Note  the.  stresses  in  the  1/2  kbar  case  are  roughly  twice 
those  in  the  1/4  kbar  case;  both  figures  reflect  a  higher  c 
for  y  =0.1  because  y  =  1.0  allows  more  plastic  flow  and 


thus  lower  stress.  This  larger  amount  of  plastic  flow  (y  =  1.0) 
is  also  responsible  for  a  larger  final  o’  at  t  =  .04/.06sec, 
because  less  elastic  strain  is  available  for  release  during 
unloading.  Figure  34  compares  load  -strain  curves  between  an 
element  with  plastic  deformation  and  one  without.  The  two 


curves  do  not  return  to  zero  load  because  P  consists  of 

max 


both  in  situ  and  weapons  effect  loading.  Note  that  the 
element  in  the  plastic  zone  is  straining  3  to  4  times  more  than 
the  elastic  element.  Although  a  large  portion  of  the  plastic 
element's  strain  is  recovered,  some  irreversible  strain  is 
present.  This  plastic  strain  is  caused  by  slippage,  crushing, 
and  rearrangement  of  rock  particles. 


Figures  35-38  compare  stresses  along  a  45°  radial. 


for  the  two  load  cases  with  varying  y  .  Figures  39-42 


compare  o  stresses  in  a  similar  fashion.  The  1/4  kbar 


stri-sses  of  Figs  35  and  36  converge  as  ::  stance  away  from  ti. 
cutout  increases,  because  only  15%  of  the  vertical  load  is 


V 


RAD  I  I 


Fig  40.  o  as  a  Function  of  Time  Along  45° Radial 
x  1/4  kbar,  y  -  1.0 


varying  with  time.  Eighty-nine  percent  of  the  1/4  kbar 
horizontal  load  varies  with  time;  respective  percentages  for 
1/2  kbar  loading  are  72%  and  95%,  explaining  the  much  greater 
variability  in  final  stresses  (Figs  37-42).  Note  in  all  cases 
that  maximum  stresses  (at  t  =  .03  for  1/4  kbar,  and  t  =  .02 
for  1/2  kbar)  near  the  cutout  exceed  the  steady  state  stresses. 
One  would  expect  the  component  stresses  to  follow  o  ,  i.e., 
stresses  at  t  =  .03  for  the  1/4  kbar  case  (t  =  .02  for  1/2 
kbar)  snould  be  a  maximum.  However,  note  in  Figs  35  and  36 
that  t  =  .03  (1/4  kbar)  stresses  are  not  consistently  the 
highest.  This  is  another  example  of  the  redistribution  of 
stresses  afforded  by  plastic  action. 

Finally,  note  that  Figs  40  and  42  show  tension  (negative 
values)  for  arithmetic  average  stresses  between  elements 

for  y  =  1.0  ,  1/4  and  1/2  kbar  loading.  Although  some 

tension  was  developed  in  the  static  analysis,  the  tension 
developed  for  the  transient  case  is  much  more  severe.  The 
1/2  kbar  case,  for  example,  developed  91  tension  elements; 
these  instantaneous  tensile  stresses  may  cause  cracks  which 
will  not  close  due  to  in  situ  loading  c.ter  the  pressure  wave 
passes,  leading  to  zones  of  fallen  rock  in  the  tunnel. 

A  Final  Note 

bishop  (Ref  1)  suggests  that  the  Drucker-Prager  yieid 
criterion  allows  excessive  plastic  flow.  Bishop  further 


e  *  «■ 
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intimates  that  this  plastic  strain  error  is  due  to  the  fact 
that  large  volumetric  strains  are  allowed  to  occur.  From 
Eg  (2.14) 


(2.14) 


where  for  direct  (volumetric)  strains 


.vp. 

a 

z£j _ 

e  .  .  =  y 

73-1 

1  +  a 

13 

_k 

2(J')  2 

_  2 

Note  that  a  =  aJ^  +  (J2>  (yield  surface) . 

In  an  effort  to  reduce  the  amount  of  plastic 
strain,  the  a  term  was  removed  from  Eq  (2.15)  so 


•  vp 

0 

’°ii 

cij =  Y 

Ik'  - 

.2(J2)l 

(2.15) 


volumetric 

that 


(3.1) 


Note  that  by  removing  a  in  this  fashion,  it  is  still  included 
in  the  o'  term  which  triggers  the  plastic  strain  rate  calcula¬ 
tion.  The  Baker,  et  al  (Ref  11)  in  situ  loading  case  was 
analyzed  using  Drucker-Prager  in  this  form  and  compared  against 
the  earlier  analysis. 

Figure  43  shows  the  size  of  the  plastic  zone,  which  was 
the  same  for  both  forms  of  the  yield  criterion.  Table  2 
compares  a  representative  stress  value  (a  )  and  plastic 
effective  strain  for  several  elements  in  the  plastic  zone, 
where  plastic  effective  strain  is  defined  as  (Ret  9) 

VP  -  fI~P  P  ~ 

’eff  J3  cijLij 
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followed  a  similar 


Table  2  shows  that  the  stresses  (o  ,  o  ,  u 

x'  z'  xy 

pattern)  as  well  as  most  of  the  effective  plastic  strains  were 


less  for  Eq  (3.1)  than  for  the  complete  Drucker-Prager  equation. 
The  validity  of  this  form  of  the  Drucker-Prager  equations  can 
only  be  proven  experimentally.  However,  the  fact  that  smaller 
effective  strains  are  predicted  is  attractive. 


TABLE  2 

Comparing  Eq  (3.1)  to  Complete  Drucker-Prager  Equations 


Element  in 
Plastic  Zone 


Eq  (3.1) 


Drucker-Prager 


1157 


1573 


1045 


1910 


1359 


1142 


1512 


. 00204 


.00183 


.00124 


.00135 


.00045 


00030 


.00183 


.00158 


.00112 


. 00083 


1328 


1785 


1185 


20  31 


1449 


1313 


112b 


00367 


00329 


00159 


00137 


00036 


00028 
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OTE :  See  Pig  43  for  Plastic  Zone 
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Summary  and  Conclusions 


A  two-dimensional  plane  stress/plane  strain  finite 
element  program  called  VISCO  was  used  to  model  the  stress 
field  around  a  deep  underground  tunnel.  The  constitutive 
equations  employed  in  VISCO  were  in  the  form  of  an  elastic- 
perfectly  plastic  material  model  developed  by  Drucker  and 
Prager.  Plasticity  was  generated  using  the  Drucker-Prager 
equations  incorporated  in  the  Malvern  flow  law.  The  following 
statements  and  conclusions  are  based  on  the  analysis  presented 
in  this  paper. 

1.  Viscoplastic  theory  can  be  used  to  model  steady 
state  elasto-plastic  material  behavior  by  setting  the  fluidity 
constant  in  the  strain  rate  equation  to  a  large  value.  Results 
obtained  in  this  effort  for  verification  purposes  agreed 
excellently  with  published  results. 

2.  It  is  possible  to  follow  the  effects  of  a  nuclear 
blast  using  viscoplastic  theory  for  both  static  and  quasi¬ 
transient  loading  of  tunnels  buried  in  tuff.  Results  produced 
for  the  steady  state  conditions  yielded  physical  activity 
recognized  in  field  tests  and  experimental  work.  The  transient 
analysis  lends  an  appreciation  of  possible  plastic  zone  growth 
allowing  for  actual  experimental  modification. 

3.  Tunnel  displacements  diffe.  -ai  by  extreme  amounts 
as  the  load  varied  in  the  quasi-static  analysis.  Spring  Line 
displacements  were  nearly  eight  times  greater,  and  crown 


k 

t. 

I. 

:« 


« 

►  ■* 


x  o» 

1 


< 


displacements  were  over  three  times  greater  for  1/2  kbar 
versus  1/4  kbar  loading.  Deformations  will  increase  as  the 
size  of  the  plastic  zone  increases. 

4.  In  the  quasi-transient  analysis,  the  size  of  the 
plastic  zone  is  a  function  of  both  loading  and  fluidity 
constant.  The  two  fluidity  constants  used  here  indicate 
trends  which  may  be  compared  against  experimental  results. 
Because  the  larger  fluidity  constant  r.  lluv:s  more  plastic 
flow,  less  elastic  strain  is  released  upon  unloading,  creating 
larger  final  (t  =  .06,  t  =  .04)  tunnel  displacements.  For 
all  cases  (both  y  and  at  any  time)  the  quasi-transient  model 
experienced  less  tunnel  deformation  than  the  comparable 
loading's  static  analysis. 

5.  Inaccurately  large  strains  can  be  developed  by  any 
elastic-perfectly  plastic  material  model.  Reports  show  that 
the  Drucker-Prager  model  may  also  predict  excessive  hydro¬ 
static  strains.  A  "modified"  Drucker-Prager  material  model 
compared  favorably  against  the  actual  equations  by  predicting 
smaller  plastic  strains. 
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